Preprint typeset in JHEP style - HYPER VERSION IFUP-TH 200; 

Topological susceptibility from the overlap 



Luigi Del Debbio 

Dipartimento di Fisica dell'Universita di Pisa and I. N.F.N. 
Via Buonarroti 2, 1-56127 Pisa, Italy 



E-mail: ldd@df . unipi . it 



Claudio Pica 

Dipartimento di Fisica dell'Universita di Pisa and I. N.F.N. 
Via Buonarroti 2, 1-56127 Pisa, Italy 



E-mail: claudio . pica@df . unipi . it 



Abstract: 

The chiral symmetry at finite lattice spacing of Ginsparg- Wilson fermionic actions 
constrains the renormalization of the lattice operators; in particular, the topological 
susceptibility does not require any renormalization, when using a fermionic estimator 
to define the topological charge. Therefore, the overlap formalism appears as an 
appealing candidate to study the continuum limit of the topological susceptibility 
while keeping the systematic errors under theoretical control. We present results 
for the SU(3) pure gauge theory using the index of the overlap Dirac operator to 
study the topology of the gauge configurations. The topological charge is obtained 
from the zero modes of the overlap and using a new algorithm for the spectral flow 
analysis. A detailed comparison with cooling techniques is presented. Particular 
care is taken in assessing the systematic errors. Relatively high statistics (500 to 
1000 independent configurations) yield an extrapolated continuum limit with errors 
that are comparable with other methods. Our current value from the overlap is 
X 1/4 = 188 ± 12±5MeV. 
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1. Introduction 

The QCD lagrangian with Nf massless quarks has a XJ(Nf)y <8>\J(Nf)A symmetry at 
the classical level. The SU^iV/)^ subgroup is spontaneously broken, yielding Nj — 1 
massless Nambu-Goldstone bosons with J p = 0~. For Nf = 3 and small quark 
masses, the pseudo NG bosons form the light pseudoscalar mesonic octet observed 
in the hadronic spectrum. If the remaining U(1) J 4 subgroup were a symmetry of the 
theory, one should observe parity doublets, while the spontaneous breaking of this 
symmetry would require the existence of a light singlet state, the rj', such that ||T|: 

m v > < VSm^. (1.1) 

The fact that both scenarios are ruled out by the spectrum of light mesons is known 
as the U(l) problem, whose solution is due to the anomalous non-conservation of the 
flavor-singlet axial current ||: 

d,jl(x) = 2N f q(x) (1.2) 
where = ■0(x)757 M ^(z), q(x) is the topological charge density: 

^) = T7T2 Tr F WHz) = d » K » (1-3) 
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and is the current: 



K ^ x ) = Sl^Ki?) (d p A a a (x) - \gr bc A b p (x)A c a (x) ) . (1.1) 



Using the large N c expansion of the theory || , where N c is the number of col- 
ors, the rf mass is related to the topological susceptibility x> through the Witten- 
Veneziano (WV) formula [[| [|: 



m v' = -jix ( L5 ) 

J 7T 

where f n is the pion decay constant, and 

X = J d 4 x (T q(x)q{0))\ YU (1.6) 
is meant to be computed in the pure Yang-Mills theory. At this stage, the definition 



in Eq. (L6) is still rather abstract, since a prescription to define the correlator as 
x — > is necessary. However, Eq. ( |1.5| ) shows already some characteristic properties: 
If x goes to a constant in the large N c limit, and using the fact that f n ~ \fW c @, one 
obtains mi, ~ l/iV c , consistently with the anomalous breaking being suppressed in 
that limit. Moreover, the WV relation connects meson spectroscopy to the nontrivial 
topological properties of Yang-Mills fields; therefore it provides a stringent test of 
the nonperturbative dynamics of QCD, which can be performed by computing the 
topological susceptibility using numerical simulations of the theory defined on a 
lattice. For the WV formula to hold, the contact term in the definition of the 
topological susceptibility should be such that 0: 

X = f d A xd,(K,{x)qm (1.7) 



A lattice topological charge density qh(x) can be defined by requiring that the 
continuum definition is recovered in the naive continuum limit, where the lattice 
spacing a — ► 0. In the quenched theory, the lattice operator is related to the contin- 
uum one through a multiplicative renormalization 0, ||: 

q L (x) = Z((3)a\(3)q(x) + 0(a 6 ). (1.8) 

where both the lattice spacing a and the renormalization constant Z are functions 
of the coupling (3 = 2N c /g 2 . The lattice topological susceptibility becomes: 

XL = X>i(aO«fc(0)> (1-9) 

X 

however, in order to compare with Eq. ( |L5| ), the divergence for x —>■ has to be 
properly renormalized, which in principle involves an additive renormalization to 
take into account the mixing with lower dimensional operators ||: 

XL = Z 2 (P)a\(3) X + M((3) (1.10) 
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Different methods have been devised to deal with the renormalization of the topo- 
logical susceptibility. Precise numerical results have been obtained by computing 
numerically Z and M non-perturbatively |[10| , IT| , |T^j or by using a cooling tech- 



nique O, QJ, 15, 16 1 . In the first case, assuming that the UV modes thermalize more 



quickly than the long-range topological ones, the renormalizations are determined 
by performing a few Monte Carlo updates on given semiclassical configurations. In- 
deed, the topological charge autocorrelation time tq is known to increase faster than 
the autocorrelation time for the Gaussian high-momentum modes as the continuum 
limit is approached, which should guarantee that the short-distance contributions 
are thermalized without altering the topological content of the lattice gauge field 
configuration. Such a behaviour has been confirmed by recent simulations, showing 
that tq actually grows exponentially with the correlation length [|16j. Similarly, the 
cooling technique also relies on the hypothesis that the topological and UV degrees of 
freedom decouple: in this case one wants to "switch off" the high-momentum modes, 
thereby obtaining Z — > 1 and M — > 0, without modifying the topological content of 
the theory. The method relies on the hypothesis that, close to the continuum limit, 
the topology of a gauge configuration is determined by the structure of the fields 
on length scales p larger than the cut-off a. If this were indeed the case, one could 
"smooth" the gauge configurations on distances d such that a<dC/), e.g. by local 
minimizations of the action. Along the cooling process two different phenomena can 
occur. Instanton-antiinstanton annihilations can take place, leading to a modifica- 
tion of the topological charge density, but harmless for the total topological charge. 
Besides, a shrinking of the instanton is also observed, which leads to a loss of topo- 
logical charge. Metastable results are obtained if indeed a <ti p. Despite the success 
of the cooling method in its various versions, it still introduces a source of systematic 
error in the measurements, which is expected to be larger on coarser lattice, or at 
finite temperature. 

Lattice formulations of fermions with Dirac operators that satisfy the Ginsparg- 
Wilson (GW) relation O] have been introduced in recent years [18|, [L^, p0| . They all 



have a global chiral symmetry [£1]] and satisfy an index theorem |22[ at finite cutoff, 
which allows the topological charge of a field configuration to be computed from the 
number of their zero-modes. In particular, the WV formula has been derived on 
the lattice: using GW fermions, one obtains a definition of the lattice topological 
susceptibility to use in Eq. ( |1.5| ) which does not need any subtraction J33J. More 
details about the computation of the topological susceptibility from the index of the 
overlap Dirac operator are presented in Sect. ||. 

Studies of the topological structure of the QCD vacuum using a fermionic def- 
inition of the topological charge have already been performed both at zero and at 
finite temperature p4|, |25|, |27|, |28|, |30], |31fl . In this paper, we investigate the 
continuum limit of the topological susceptibility in the overlap formulation of lat- 
tice QCD; using high statistics and a detailed comparison with the cooling method, 
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we assess some of the systematic errors that appear in the computation and in the 
extrapolation. The results of our simulations are presented and discussed in Sect. |[ 
Finally, we summarize our results and discuss some outlooks in Sect. f|. 

2. Lattice formulation 

2.1 Witten-Veneziano formula 

The massless overlap Dirac operator can be written as: 

D = -[l + l5 e(H w (m))} (2.1) 
a 

where H^(m) is the Hermitian Wilson-Dirac operator 7 5 Z>vv( — m ) with mass pa- 
rameter —77i, e is the sign function, and a = a/m. It is easy to check from Eq. ( |2.1| ) 
that D verifies the GW relation: 



7 5 L> + L>7 5 = aDj 5 D, 



[2.2) 



In principle, m can be chosen at will in the range < m < 2, with the different 
overlap operators defined in this way yielding the same physics in the continuum 
limit. In practice, both the spectrum of H^(m) and the locality of the overlap 



depend on m p2\ \. The dependence on m of the topological susceptibility is studied 
in detail and presented below. 



In deriving the WV relation on the lattice using GW fermions, one obtains [23 



f 2 

J TV 

2N-"" r > 



2 -4 

m„, = a 



(^Tt[ 15 D(x,x)] ^Tr [75^(0,0)]) 



a 

2a 



YM 



(2.3) 



A comparison with Eqs ( |1.5|) and (|1.6| ) yields for the lattice topological density: 



<Il{x) = — Tr [y s D(x,x)] 



Hence, for the topological charge: 



n_ — n_ 



(2.4) 



(2.5) 



where n± are the numbers of zero modes of the overlap Dirac operator with positive 
and negative chirality respectively. Using Eq. ( |2.5| ), we get for the susceptibility on 
a finite lattice of volume V: 



XL 



((n_-n + ) 2 ) 
V 



(2.6) 



and xl = ^ 4 X) without additive nor multiplicative renormalizations needed. 
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2.2 Index of the overlap 

Let us now briefly summarize the main features of our implementation of the overlap 
operator. In order to study topology, we only need to concentrate on two aspects: the 
application of the overlap operator to a generic spinor field ip, and the computation 
of the index of D, which requires the study of its zero-modes. Equation ( p . 1|) can be 
rewritten as: 



a 



m 



(2.7) 



so that the implementation of the sign function e{H^{m)) is made explicit. The 
Neuberger operator D is invariant under a rescaling of the Hermitian Wilson-Dirac 
operator H^(m); the normalization used in this paper is such that the spectrum of 
Hw(m) is bounded by 1. The implementation of the inverse square root of a sparse 
matrix with a potentially large condition number is a demanding task and the interest 
in the overlap formulation of lattice QCD has triggered several studies devoted to 
its implementation, and a "standard" technique has been established. Some of the 
lowest-lying eigenvalues of H w (m) are computed exactly, and the inverse square root 
is trivially defined in the basis provided by the corresponding eigenvectors. In the 



orthogonal subspace, 1/ a/ H^(m) is computed either using a polynomial expansion 
or a rational approximation. Detailed discussions and an exhaustive list of references 
can be found in Refs. |33|, 

In this work, the 15 lowest eigenvalues of H^(m) are computed using an accel- 
erated conjugate gradient algorithm []35| and the corresponding low- lying modes are 
treated exactly. The eigenvalues of ifj^(m) are required to have an absolute preci- 
sion of at least 10" 15 . The conjugate gradient search stops when all the eigenvalues 
have the desired precision. The error on the eigenvalues is estimated exploiting the 
quadratic convergence of the conjugate gradient algorithm in a way similar to that 
of Ref. P5| . Since the quadratic convergence regime is only reached asimptotically, 
we used the rigorous bound \g\ as an estimate of the error as long as it stays bigger 
than 10~ 4 . Then the error is assumed to be 

K^4, (2.8) 
1 - r ' 11 



where primed variables refers to previous step values, /x is the Ritz functional value 
and g its gradient. The numerical constant factor K was fixed from preliminary 
study comparing the true error and the error estimate. The numerical value of K is 
10. The a priori bound \g\ results to be about a factor 1000 greater than the above 
estimate. 

The approximation on the orthogonal subspace is done using Chebyshev poly- 
nomials, so that the inverse square root reads: 

l/y/H^(m) = J2 -TxtMM + E ^T n {H^{m)) (2.9) 
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where Aj, V{ are the eigenvalues and eigenvectors that have been computed exactly, 
and T n are Chebyshev polynomials. Equation ( |2.9|) summarizes the algorithm that 
evaluates e(H w (m))il), and hence Dip. The gauge covariance, 75-hermiticity, and 



locality of our overlap operator have been succesfully tested ||32 . 

In order to compute the index of the overlap operator, one has to compute 
the eigenvalues of D^D. Since [75, D^D^ = 0, one can work in subspaces of given 
chirality. Let us denote by P± the projectors on the chirality subspaces, using the 
GW relation one can show: 

D ] DP ± = P ± D ] DP ± = -D ± (2.10) 

a 

where D± = P±DP±. Hence, by working in fixed chirality subspaces, the compu- 
tation of (D^D)ip only requires one application of the sign function. Moreover, one 
can show that the spectra of D± exactly coincide, except for the number of zero 
modes in the two sectors. Actually, configurations with zero modes in both sectors 
are statistically irrelevant; excluding such a possibility simplifies the numerical task. 



Following the strategy outlined in |34j], we run the minimization program in both 
sectors simultaneously, until one of them is identified as a sector in which there are 
no zero modes. We require the eigenvalue /1 to satisfy the bound \i — \g\ > 0, and 
the relative precision of the eigenvalue to be at least 5%. A refined search is then 
performed in the other chirality sector only, in order to count the zero modes. A 
standard conjugate gradient algorithm is used to find the smallest eigenvalue. This 
procedure is repeated until an eigenvalue compatible with that in the other sector 
is found. We use this value as an upper bound for finding the zero modes. A zero 
mode is identified with an eigenvalue incompatible within its error with this upper 
bound and whose magnitude is less than 20% of the upper bound. The search for 
zero mode ends when the current eigenvalue differs from the upper bound by less 
than 10% and its relative precision is at least 10%. 

2.3 Spectral flow 

If D is defined using Hw(m) for some given value of m, then one can show that 
its index corresponds to the number of eigenvalues of Hyf(fj) that cross zero (level 
crossings), for < /J < m JT9|]. The spectrum of H^(fi) is characterized by three 



regions: for small //, < /j < mi, the Wilson-Dirac operator describes fermions with 
positive physical mass and hence the spectrum does not show any level crossing; as 
/j is increased a second region mi < fi < m 2 exists, where the spectrum is gapless 



and crossings occur; finally the gap should open again for \i > m 2 fl26fl . One would 
expect that mi,m 2 — > 0, as the continuum limit is approached. For the values of (3 
typically used for numerical simulations, the gap closes at some value of /1 ~ 0.8 and 



remains closed until the end of the allowed region /x = 2. It was pointed out in 
that the zero modes that appear far from mi have a size of the order of the lattice 
spacing and should not affect the physics in the continuum limit. 
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To determine the number of level crossing in the spectrum of the hermitian 
Wilson-Dirac operator as a function of the mass parameter m we performed a spectral 
flow analysis. Using the information about the spectrum of H w (m) f3~?fl , we introduce 
a new procedure to locate the crossings in a given mass region [m a ,m u ]. Starting 
with m = m a , we find the two smallest (in modulus) eigenvalues of H^(m) together 
with their derivatives with respect to m; the value of the mass m is then increased 
by 5m, and two new eigenvalues and derivatives are computed. The algorithm to 
identify a level crossing in the interval [m', m] (primed variables indicate previous 
step values), is the following: 

1. find the two smallest eigenvalues (Ao < Ai) and eigenvectors of H^{m); 

2. calculate and store eigenvalues derivatives, together with respective eigenvalues; 

3. find which of the two eigenvalues computed in step (1), Ao or Ai, can most 
probably be considered the continuation of A . This is done by comparing the 
average of the derivatives with the incremental ratio of the eigenvalues; 

4. a level crossing between [mf, m] is identified if the eigenvalue defined in step (3) 
has a different sign from A and the eigenvalue derivatives match within 10 deg; 
the exact crossing value and derivative are estimated by linear interpolation; 

5. increment the mass by 5m — \X\\ and if m < m u repeat from step (1) 

The eigenvalues at step (1) are found using the same accelerated conjugate gradient 
algorithm cited before with a numerical precision of 10 -12 . If at point (4) there is 
the possibility of a level crossing but the derivatives do not satisfy the desidered 
constraint we tag the interval and repeat the procedure on that interval with a 
smaller mass step size. The value of 5m in step (5) is chosen so that only the 
lowest eigenvalue should possibly have crossed before m', as granted by the bound 



on derivatives \d\/dm\ < 1 [37 . 



3. Numerical results 



3.1 MC data 



In order to study the continuum limit of the topological susceptibility, we performed 
simulations of the SU(3) lattice gauge theory in the standard Wilson formulation, 
at three values of the coupling (3 = 5.9, 6.0, 6.1. The lattice volume are 12 4 , 12 4 , 16 4 
respectively, so that each lattice has a linear size L > 1 fm. The lattice spacing 
can be set either from the string tension or from the low-energy reference scale 
ro [36|. The relevant numbers for our simulations are summarized in Tab. |], where 



the dimensionless ratio r y/a is also reported. The latter only varies by a few percent 
over the range of (3 that we explore, indicating that the two scales yield compatible 
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p 




r /a 




a (T (fm) 


a ro (fm) 


5.9 


0.2605(14) 


4.48 


1.16704 


0.12 


0.11 


6.0 


0.2197(12) 


5.37 


1.17979 


0.10 


0.09 


6.1 


0.1876(12) 


6.32 


1.18563 


0.09 


0.08 



Table 1: Lattice spacing in units of the string tension and of the reference scale ro for 
the values of (3 used in this work. The physical values a a and computed assuming 

\fo = 440 MeV and ro = 0.5 fm respectively. 



results in the continuum limit. Both quantities will be used to build scaling ratios 
and study the continuum limit of xl- 

For each value of (3, we compute the topological susceptibility using the index 
of the overlap, and we study its dependence on the parameter m by a spectral flow 
analysis. Furthermore, we also use a cooling algorithm to evaluate the topological 
charge on the same configurations and compare the outcomes. Since the lattice 
artifacts in the two definitions can be completely different, there is no reason a priori 
for these two quantities to agree at finite lattice spacing. However, one would like to 
check that the two methods do agree as the continuum limit is approached. Moreover, 
the "systematic errors" of the cooling technique being relatively well understood, such 
a comparison, even at finite values of (3, can shed some light on the effectiveness of 
the fermionic method in detecting topology. 

The gauge configurations are generated by a Cabibbo-Marinari algorithm, alter- 
nating heatbath and microcanonical updates in a ratio of 1:4. In order to guaran- 
tee the statistical independence of our configurations, measurements are separated 
by 1000 updates, which is much larger than the estimated autocorrelation time, 
tq < 200, at these values of the coupling [TI|. On each configuration the topo- 
logical charge is computed using cooling, a spectral flow analysis, and the index of 
the overlap. The number of configurations analyzed for each lattice is reported in 
Tab. |2|, together with the values of the topological charge and susceptibility. It is 
clear from the raw lattice data, that the average topological charge is always zero 
within errors, as it should be. Moreover, as we shall see in more detail in what 
follows, the differences between the different methods do decrease as the continuum 
limit is approached. The simulations have been performed on a cluster of Pentium-4 
Xeon processors at 2.2 GHz, using SSE2 instructions to implement the most time- 



consuming operations [3~8"|. The total time used for the actual simulations is about 
25 CPU-months; with a breakdown of 15 CPU-months for the overlap and 10 for the 
spectral flow. 

The time-history of the topological charge for several values of (3 and m are 
shown in Fig. p], [|, where one can see that the algorithm does tunnel from one 
topological sector to the other and that the distribution of values is symmetric and 
peaked around zero. It is also clear from the figures that the separation between the 
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lattice 


conf 


method 


m 


Q 


IQI 


X (xlO 4 ) 


5.9 


12 4 


431 


overlap 


1.0 


0.05(7) 


1.21(5) 


1.27(9) 






510 


overlap 


1.5 


0.17(9) 


1.51(5) 


1.83(11) 






432 


cooling 




0.01(8) 


1.34(5) 


1.43(9) 








Ref. 1161 
i_ i 








1.544(7) 


6.0 


12 4 


956 


overlap 


1.0 


0.03(4) 


0.90(3) 


0.70(3) 






956 


overlap 


1.5 


0.01(4) 


0.98(3) 


0.80(5) 






956 


cooling 




0.03(4) 


0.97(4) 


0.72(3) 








Ref. pi 








0.728(5) 


6.1 


16 4 


440 


overlap 


1.0 


0.09(8) 


1.12(4) 


0.33(2) 






440 


overlap 


1.5 


0.09(8) 


1.15(5) 


0.35(3) 






355 


cooling 




0.09(8) 


1.13(5) 


0.34(5) 








Ref. [16 








0.382(6) 



Table 2: Topological charge Q and susceptibility \l- 




Figure 1: The time history of Q for (3 = 5.9, m=l (bottom), and m=1.5 (top). 

measurements is large enough to provide statistically independent data. At larger 
(3, the topological charge distributions are narrower, as one would expect from semi- 
classical arguments. Comparing the histograms at different values of m, we see that, 
as m is increased, the charge distribution broadens. 

The mass dependence of the topological susceptibility, obtained from a spectral 
flow analysis, is shown in Fig. ||. We observe a large effect at lower values of /3, where 
the value of topological susceptibility as a function of m shows a plateau within the 
statistical errors only for m <; 1.5. For m ~ 1.0, the value of the susceptibility is 
about 75% of the asymptotic result. For /3 > 6.0, the same plateau starts for m < 1.0. 
Such dependence should be taken into account in studies performed at fixed value of 
m. The density of crossings vs. m is reported in Fig. Q the expected behaviour is 
clearly observed. The value of mi, where the gap closes, decreases as the continuum 
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Figure 2: The time history of Q for (3 = 6.0, m=l (bottom), and m=1.5 (top). 
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Figure 3: Topological susceptibility xl as a function of m for (3 = 5.9 (top), 6.0 (middle), 
and 6.1 (bottom). 



limit is approached; in the same limit, the level crossings cluster in the vicinity of mi, 
so that the topological susceptibility becomes effectively m-independent. For lattice 
sizes L > 0.9fm, our results compare well with those presented in Ref. p6[ , where 
finite size effects were observed on smaller lattices, e.g. on a 8 3 x 16 at (3 = 6.0. 

The data obtained from the overlap can be compared with those obtained by a 
cooling algorithm on the same configurations. Data in Tab. |2| show that the discrep- 
ancy between the fermionic and bosonic estimator of the topological charge (Q / and 
Q g , respectively) decreases as the continuum limit is approached. The distribution 
of AQ = Q g — Qf is reported in Fig. || for m = 1 and (3 = 5.9 and 6.0. For (3 = 6.0, 
the agreement of the two methods at m = 1 is around 90%; a comparison of the two 
values of (3 shows that the discrepancy between the two methods decreases as the 
continuum limit is approached, consistently with the scenario outlined in Ref. pi]]. 
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Figure 4: Distribution of level crossings as a function of m for (3 = 5.9 (bottom), 6.0 
(middle), and 6.1 (top). 



The agreement between the two results is also reassuring as far as finite size effects 
are concerned. Indeed, the high statistics results obtained in Ref. flTfJ with a cooling 



technique show that the finite size effects on the topological susceptibility at the 
values of (3 considered in this work are of the order of a few percent. If the agree- 
ment between the two methods is not a coincidence, the cooling and the fermionic 
estimators are expected to display similar finite size effects, which are known to be 
well below the statistical accuracy achieved in this work. 

In Fig. we report the rate of discrepancies between the fermionic and gluonic 
definitions as a function of m, for the three values of the coupling. It is interesting to 
note that the best agreement is obtained at intermediate values of m. The level cross- 
ings in the vicinity of the threshold m\ are the ones that build the topological charge 
detected by the cooling method. Level crossings further away from m\ correspond 
to topological excitations that are eliminated by the cooling process, suggesting that 
they are associated to modes at the scale of the cut-off. 

3.2 Extrapolation to the continuum limit 

In order to study the extrapolation to the continuum limit, we build the adimensional 
ratios C a = Xl/& 2 and C ro = XL r o, and report their values in Tab. [3[ To compare 
directly with other works, we also compute the fourth root of the scaling ratios. The 
scaling ratio C a is displayed in Fig. ^ as a function of a 2 a, since the leading scaling 
corrections are expected to be 0(a 2 ). The data at m = 1 do not display any sizeable 
dependence on a 2 , and a naive extrapolation can be obtained simply by quoting an 
interval such that all data points are included: 

C a (p -f oo) = 0.029 ± 0.05 (3.1) 
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Figure 5: Distribution of the difference AQ = Q g — Qj for m = 1, ft = 5.9 (left), and 6.0 
(right). 



Consistently with the fact that the m-dependence becomes negligible as ft is in- 
creased, as shown in Fig. |3|, the data point at m = 1.5 and ft = 6.1 is included in the 
above interval. 

Alternatively, one can take a continuum limit by fitting the data either to a 
constant, or a linear function. We obtain: 

C a (ft -> oo) = 0.028 ± 0.02 (constant) (3.2) 
C a (ft -f oo) = 0.026 ± 0.04 (linear) (3.3) 

(3.4) 

where the errors on the fitted parameters are estimated using a bootstrap method. 
The results obtained are in excellent agreement with previous determinations. 

It is clear from Fig. |^that there are sizeable scaling violations for m = 1.5, such 
that the extrapolated results in the continuum limit are dominated by the large value 
at ft — 5.9. In order to get a more reliable result one should increase the statistics 
for the current data and add another point closer to the continuum limit, e.g. at 
ft = 6.2. This would require a larger lattice in order to avoid finite size effects, an 
effort that is beyond the scope of this work. 

Data obtained from overlap operators with different values of m can be put 
together in a single constrained fit, assuming that they all yield the same continuum 
limit, and allowing for quartic terms in the extrapolation: 

C a , m {ft) = C a {{3 - oo) + c ltm a 2 (ft) + c %m a\ft) (3.5) 
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Figure 6: Discrepancy rate between the gluonic and fermionic determinations of the 
topological charge as a function of m for j3 = 5.9 (top), 6.0 (middle), and 6.1 (bottom). 

with Co independent on m. Such a procedure allows a non-linear dependence on a 2 
to be taken into account; however, one should remember that the data used for the 
simultaneous fit are not independent, so that the bootstrap errors on the fitted values 
should only be taken as an indication. The result of the constrained fit is: 



C a {(3 



oo 



0.025^ 



(3.6) 



As a consistency check, one can fix C a (/3 — > oo) to the value obtained above, e.g. 
by fitting the data at m = 1 to a constant, and fit the data for both values of m to 
Eq. (|3.5|) with only four parameters. The data fit this ansatz with a x 2 /dof ~ 1. 



The same analysis can be performed for the scaling ratio C ro . The reference 
scale tq has been computed to great accuracy and a formula for the interpolation 
to arbitrary values of j3 G [5.7,6.57] does exist J36]. The use of r for the scaling 
analysis guarantees a greater uniformity when trying to unify data from different 
studies The results of the fits are: 



C ro (P -> oo) = 0.054(3) for m = 1.0, constant fit 
C ro (P -> oo) = 0.055(8) for m = 1.0, linear fit 



(3.7) 
(3.8) 



A naive extrapolation, which includes all the results for m — 1, yields C ro 
0.055(10). 

Combining the two determinations, we obtain: 



X 



1/4 = 188 ± 12 ± 5 MeV 



(3.9) 
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method 


m 


c a 


cT 


' u 


Wl/4 


5.9 


overlap 


1 


0.0273(23) 


0.407(9) 


0.0509(32) 


0.475(8) 




overlap 


1.5 


0.0390(37) 


0.445(10) 


0.0727(52) 


0.519(8) 




cooling 




0.0308 26 


0.419(9) 


0.0574(36 


0.489 8 




Ref. 1161 
i_ i 




0.0334(9) 


0.428(3) 


0.0622(4) 


0.499(1) 


6.0 


overlap 


1 


0.0300(19) 


0.416(7) 


0.0581(25) 


0.491(5) 




overlap 


1.5 


0.0343(25) 


0.430(8) 


0.0664(33) 


0.508(6) 




cooling 




0.0311(19) 


0.420(7) 


0.0602(25) 


0.495(5) 




Ref. |16| 




0.0312(7) 


0.420(3) 


0.0604(4) 


0.496(1) 


6.1 


overlap 


1 


0.0260(26) 


0.402(10) 


0.0515(38) 


0.476(9) 




overlap 


1.5 


0.0280(28) 


0.409(10) 


0.0555(42) 


0.485(9) 




cooling 




0.0282(27) 


0.410(10) 


0.0558(41) 


0.486(8) 




Ref. 16 




0.0308(7) 


0.419(4) 


0.0611(10) 


0.497(2) 



Table 3: The scaling ratios C a , Cj , C ro , and C ri 



where the first error is statistical and the second corresponds to the two different 
methods of setting the scale. 

4. Conclusions 

Using a fermionic estimator to compute the topological susceptibility in the contin- 
uum limit has some appealing features: the chiral symmetry at finite lattice spacing 
fixes the renormalization of the lattice operators, so that no renormalization constants 
are needed. This is also the case for determinations based on the cooling technique; 
nonetheless the latter are based on assumptions on the separation of gaussian and 
topological modes. 

However, the fermionic method also introduces systematic errors, which one 
would like to control: the dependence on the mass parameter m, that appears in 
the definition of the overlap, finite size effects, and the size of the lattice artifacts at 
current values of the coupling. 

In agreement with theoretical arguments and previous numerical investigations, 
we find that the topological suceptibility does vary with m at small values of /3, 
and that this dependence disappears as the continuum limit is approached. Studies 
using the overlap at fixed value of m should take this source of systematic error into 
account, even though we expect the magnitude of the systematics to depend on the 
observables. Finally, comparison with results obtained using a cooling procedure on 
the same configurations shows that the two determinations agree in the continuum 
limit. For values of m > 1.2, the fermionic method detects topological charges that 
are smoothed away by the cooling procedure, which could be related to degrees of 
freedom at the cut-off scale. 
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Figure 7: The scaling ratio C a as a function of a 2 cr. The grey band is the result of the 



continuum extrapolation in Ref. [16]. 



On our lattices of linear size L ~ lfm, the comparison with cooling suggests that 
finite size effects are small compared with our current statistical error. 

Our results for the continuum extrapolation are compatible with other determi- 
nations. In order to further constrain our extrapolation, one more point at /3 — 6.2, 
larger lattices, and higher statistics are needed. Moreover, improved actions could be 
used to reduce lattice artifacts and possibly achieve a faster convergence towards the 
continuum limit. Based on our experience, a precise determination of the continuum 
limit extrapolation is well within reach. 

If one is only interested in the topological charge and susceptibility, a spectral 
flow analysis is sufficient and can be implemented quite efficiently using the procedure 
outlined above. Instead studies of the overlap operator are needed if one wants to 
investigate the correlator of the topological density (^lO^i^O)). 
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